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ABSTRACT 


Panel Flutter and Sonic Fatigue Analysis for RLV 

Chuh Mei and Guangfeng Cheng 
Department of Aerospace Engineering 
Old Dominion University, Norfolk, Virginia 23529-0247 

A methodology is presented for the flutter analysis of the seal of thermal protection 
system (TPS) panel of X-33 Advanced Technology Demonstrator test vehicle. The seal is 
simulated as a two-dimensional cantilevered panel with an elastic stopper, which is 
modeled as an equivalent spring. This cantilever beam-spring model under the 
aerodynamic pressure at supersonic speeds turns out to be an impact nonlinear dynamic 
system. The flutter analysis of the seal is thus carried out using time domain numerical 
simulation with a displacement stability criterion. The flutter boundary of the seal is 
further verified with a family of three traditional and one nontraditional panel flutter 
models. The frequency domain method that applies eigenanalysis on the traditional panel 
flutter problem was used. The results showed that the critical dynamic pressure could be 
more than doubled with properly chosen material for the base stopper. The proposed 
methodology can be easily extended to three-dimensional panel seals with flow 

angularity. 



FOREWORD 


This technical report contains the research results on supersonic panel flutter analysis of 
thermal protection system (TPS) seals. The research is sponsored by a NASA Langley 
Research Center Grant NAG- 1-2 150, entitled '‘Panel Flutter and Sonic Fatigue Analysis 
for RLV.” Technical monitor is Dr. Max L. Blosser, Metal and Composite Thermal 
Structures Branch, NASA Langley Research Center. Two journal publications have 
resulted from this research grant: (1 /‘Nonlinear Response of Composite Panels Under 
Combined Acoustic Excitation and Aerodynamic Pressure,” A LA A Journal, Vol. 38, No. 
9, 2000, pp. 1534- 1542. and (2)“A Methodology for Supersonic Panel Flutter Analysis of 
TPS Seals,” accepted for publication and to appear in Journal of Aircraft, July 2001. The 
authors are extremely thankful to Max Blosser for many technical discussions and 
suggestions. 
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NOMENCLATURE 
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= element and system aerodynamic influence 
matrices 
= beam width 

= aerodynamic damping coefficient 
= element nodal displacement vector 
= Young's modulus 
= impact force 

= non-dimensional aerodynamic damping 
= element and system aerodynamic damping 
matrices 

= beam thickness 
= area moment of inertia 
= equivalent spring stiffness 
= element and system stiffness matrices 
= beam length 
= Mach number 

= element and system mass matrices 
= beam element shape functions 
= aerodynamic pressure 
= dynamic pressure 

= element and system load vectors 
= airflow speed 
= transverse displacement 

= transverse deflection at beam free end 
= system nodal displacement vector 
= panel damping rate 
= non-dimensional eigenvalue 
= non-dimensional dynamic pressure 
= mass ratio 

= beam and air mass densities 
= eigenvector 

= panel and reference frequencies 
= complex panel motion parameter 


IV 
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INTRODUCTION 

Various methods have been employed to solve the panel flutter 1 " problem, which is 
among the general considerations for designing structures under supersonic/hypersonic 
conditions. In point of dynamics, the flutter boundary generally takes the form of critical 
non-dimensional dynamic pressure. The solution means can be classified in two 
categories: 3 frequency domain methods and time domain methods. While using finite 
element analysis, frequency domain methods are feasible while the problem can be cast 
as an eigenvalue problem. However, when it’s mathematically difficult to execute 
eigenanalysis, time domain methods can be applied for determining flutter. The 
acceptability of the results depends on the numerical stability of the methods adopted. 

The problem cited in this paper arises from the stability analysis of the seal of 
superalloy honeycomb thermal protection system (TPS) panel of X-33 vehicle. From the 
blueprint of the seal for two adjacent TPS panels, a sketch was made and shown in Fig. 1. 
The vibration of the covering seal is restrained with the neighboring panel. As the 
overhead airflow speed increased, flutter will occur. The seal is restrained to move 
downward by the top surface of the adjacent panel and a tiny gap (or zero gap) exists 
between the seal and the adjacent panel. Therefore, a simplified two-dimension model for 
the seal (flow angle A = 0°) is idealized as a cantilever beam with a stopper below its free 
end, as shown in Fig. 2. Accordingly, the problem turns out to be a vibration problem 
with an impact at some amplitude. 

If the stopper is very stiff so that it can be taken as a rigid stopper, the beam will 
undergo impact force during its vibration and a restitution coefficient can be included to 
allow for energy loss. 5 However, the stopper should be treated as elastic, then more 
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considerations are needed. For this case, an equivalent spring can be used to replace the 
stopper. In investigation of the family of beams with left end clamped, right end using 
various supports as shown in Fig. 3, flutter analysis of the spring stopper model as well as 
the bar stopper model, which is closer to real flutter problem, belongs to nontraditional 
panel flutter analysis. In the present paper, time domain method is developed and 
employed to solve this nontraditional flutter problem. Apparently, flutter of model I, 
model IV and model V in Fig, 3 is traditional panel flutter since eigenanalysis is 
applicable. 

The systematic procedures of finite element analysis were followed, considering the 
aerodynamic effect with quasi-steady first-order piston theory. For the spring stopper 
model, the equivalent spring stiffness is from the assumption that the stopper undergoes 
elastic deformation. Thus, the impact force is the restitution force of the equivalent 
spring. Main interest is focused on the variation of critical dynamic pressure versus 
stopper stiffness. The flutter boundary is compared with those of similar cases with 
simpler boundary conditions. 


FINITE ELEMENT FORMULATION 

The classical Euler-Bemoulli beam element can be extended with consideration of 
aerodynamic effects: additional aerodynamic damping and influence matrices need to be 
derived. The equations of motion (EOM) which govern the motion of the beam are 
derived from the principle of virtual work. 
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Element Matrices and System Equations of Motion 

From the quasi-steady First-order piston theory, the aerodynamic pressure applied on 

the beam is 
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where p = p (i V 2 / 2 is the dynamic pressure, V is the airflow speed, p a is the air density, 

is the Mach number, and E, I, b, L are property parameters of the beam. The non- 
dimensional dynamic pressure and aerodynamic damping coefficient are given by 
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where p is the mass density of the beam, coo = (EI/pAL ) " is the reference frequency and 
C a is aerodynamic damping coefficient (C a ~ p/Moo for Mco»l, and p = p a L/ph is the 
mass ratio). 

For Euler-Bemoulli element 
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The strain energy of beam element is 
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^ 0 

where h e is element length. Using Eqs. (6), the variation of strain energy can be expressed 
as 
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Similarly the virtual work done by external forces is 
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Use the principle of virtual work <5U e = 6W e 
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Thus, the element EOM is 


A- [m]{d}+ — [g ]{d}+ (A[a a ] + [K e ]){d}= {Q, } (10) 

oj q 

where the element matrices and load vector are 

M = [$ ] = 77 J L N _T L N Jd x 
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W=£/jL N I..L N i„ d * 
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After assembling, the system EOM becomes 

-4 [M ]{w}+ ^ [g]{w} + {A[A a ]+ [K 1){W}= {Q} (11) 
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Solution Procedures 

While using frequency domain analysis to determine the critical dynamic pressure tor 
traditional self-excited panel flutter problem, the homogenous form of Eq. (1 1) is used: 

± (A , ][w}+ i*. [G}{vV} + U[A,]+[^W'V}= { 0 } (12) 

0) o" CO 0 

The flutter boundary can be easily turned to be an eigenvalue problem by assuming the 
response as: 

{W} = {0}e n ' (13) 

where Q = (X + ico is generally a complex eigenvalue, or and u> are the panel damping rate 
and frequency, respectively. By substitution of Eq. (13), Eq. (12) leads to a 
nondimensional eigenvalue problem of the form 

tfM]{OJ = (A[AJ + [*]){<&} (14) 


where the nondimensional eigenvalue is 
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This is because that the mass matrix [M] equals to the aerodynamic damping matrix [G] 


as shown in Eq. (10). While there is no airflow, A = g a - 0, the eigenvalues k are real and 
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positive corresponding to the square of natural frequencies of free vibration of traditional 
panel systems (Model I, IV or V). As A is increased in value monotonically from zero, the 
symmetric, real stiffness matrix [K] is then perturbed by the skewed aerodynamic 
influence matrix [A a ] so that two of the eigenvalues approach each other until they 
coalesce to the value k„ at A = A cr . 

The aforementioned impact problem is not applicable to frequency domain analysis due 
to the difficulty of describing the boundary condition at the free end (or the impact force). 
Therefore, simulation of beam time response under airflow and impact force is hopeful to 
unveil the flutter characteristic of the cantilever beam. The equivalent spring has a 
stiffness of k $ = E S A S /L S with E s , A s , L s to be the Young’s modulus, cross section area and 
length of the stopper, respectively. Among those popular time integration schemes, the 
Newmark-/? method was adopted here for its stability. The initial conditions used are the 
first bending mode normalized with W max /h = 0.1 and zero initial velocity. Time step 
length At was decided on basis of a bisection strategy: the applied time step size is the 
one giving slightly different critical dynamic pressure compared to the Acr attained by its 
double. The ratio W max /h was checked at each time step for judgement of flutter 
boundary. The upper limit of W max /h is 4.0, i.e., simulation will be ceased as W max /h 
equals to or is greater than 4.0. The impact force is 

f 0 w,>- Gap m 

f \~ k s (|W, | “ Gap ) W e < -Gap 

where W e is the transverse deflection at the free end of the beam, Gap is the size of the 
gap marked up in Fig. 2. 
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RESULTS AND DISCUSSIONS 


To validate the finite element model is to compare the results with available analytical 
analysis results. However, it’s found that for the five members in Fig. 3, no analytical 
flutter boundaries are presented in literature. The alternative way adopted here is to check 
the flutter boundaries of clamped-clamped beam, and both ends simply supported beam 
with corresponding analytical results 6 Figure 4 shows that for these two types ot 
supporting conditions (traditional flutter analysis), finite element analysis agrees well 
with analytical analysis. 

Then, the flutter behavior of a steel beam with dimension L = 1 in.( 0.0254 m ), b = 
0.1/>z.( 0.00254 m ) and h = 0.006 in. ( 1.524x1 0' 4 m ) was investigated. The beam was 
discretized with 10 elements. To simulate the tiny gap between the beam free end and 
stopper top, a gap size of Gap/h = 0.0001 was employed. Computation was carried out for 
several distributed stiffness values within the range k s = 0-1.2x10 6 Ibf/in. ( 2.1xl0 8 N/m 
), with the upper limit corresponding to steel beam versus steel stopper case. An 
aerodynamic damping of C a = 0.02 was assumed and used in the analysis. 

Figure 5 is the critical dynamic pressure versus spring stiffness curve of the equivalent 
spring model II. It shows the existence of an intermediate stiffness ranges roughly from 
40 lbf/in.(7005 N/m) to 40,000 lbf/in. (7.0xl0 6 N/m), within which the critical dynamic 
pressures are over 300. The right and left sides adjacent 

to this range are low flutter boundary area having nearly same critical dynamic pressure 
level as that of the cantilever beam (Fig.6). The time response, force history of 
representative points A(k s =0.4, A=140), B (k s =2,000, A=300), C(k s =2,000, A=323), 
D(k s =200,000, A= 138 ), E(k s =0.4, A=160 ), F( k s =200,000, A=150 ) are shown in Fig. 
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7— Fig. li, respectively. Stable time response has a decaying or limit cycle oscillation 
response history (Fig. 7(a), Fig. 8(a) and Fig. 10(a)). Oppositely, unstable vibration 
(flutter) is featured as either gradually divergence (Fig. 9(a)) or abrupt increase of W max /h 
till the upper limit is reached soon (Fig. 1 1). 

Point A is below the flutter boundary for low spring stiffness area. Therefore, as shown 
in Fig. 7(a), the time response is stable. The spring is so soft that the stopper has little 
effect on the vibrating beam. Thus the critical dynamic pressure is very near that of the 
cantilever beam. This can also be verified by the force shape (Fig. 7(d)) of one whole 
impact. The dominating frequency standing out in Fig.7 ( b ) is obviously close to the 
coalescence frequency of cantilever beam as can be identified from Fig. 6 approximately. 
While the stiffness increased to the intermediate stiffness area, the rebounce phenomenon 
denoted by the force shapes in both Fig. 8 (d) and Fig. 9(d) complicated the impact 
process. Both the stable motion (Point B) and unstable motion (Point C) are dominated 
by several modes, instead of the coalescence frequency, as shown in Fig. 8(b) and Fig. 
9(b). However, it’s found that the stiffer the spring stopper, the fewer the rebounce. 
Point D is very close to the flutter boundary of k s = 200,000, which is Vr = 139. 
Comparing the response Power Spectrum Density plots shown in Fig. 10(b) and Fig. 7(b), 
they have very similar dominating frequency and shape. So the critical dynamic pressure 
of point D should be close to that of low spring stiffness case. The beat phenomenon 
observed from Fig. 10(a) and ( b ) implies the lowest two modes are about to coalesce. 

Figure 3 listed a family of left end clamped beams with ascending stopper stiffness. 
Among these five types of beam, the flutter boundaries of cantilever beam (model I), 
connected spring stopper beam (model IV) and right end simply supported beam (model 
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V) can be determined by frequency domain analysis. The results are shown in Fig. 6, Fig. 
12 and Fig. 13, respectively. 

The other two models, spring stopper beam and bar stopper beam, are inapplicable to 
frequency domain analysis. The bar stopper beam model is more alike the seal of X-33 
TPS panel. The analysis of this model will help with understanding of the right-hand side 
of the curve in Fig. 5, time domain analysis of it was also executed. The 0.25 in.{ 
0.00635m ) long bar stopper has a section area of 0. 1x0. 1 in' ( 6.45 lx 10 6 m" ) and is 
modeled with 5 bar elements without damping (structural damping will effect the critical 
dynamic pressure slightly). Learning of stopper material properties turns to be necessary. 
Similarly, the Newmark-/? method was used as time integration scheme. Results from 
analyses of three types of stopper material are listed in Table 1 to show the trend of 
flutter boundary variation. Obviously, the drop of critical dynamic pressure as stiffness 
increased to very stiff is confirmed. 

Table 1 Critical dynamic pressure of bar stopper beam with the stopper 

modeled as bar elements 


Stopper Material 

Boron 

Aluminum 

Steel 

Density( lb/in 3 ) 

0.083 

0.098 

0.283 

Young’s Modulus 

E s ( msi ) 

0.0638 

8.99 

30.2 

Equivalent Spring Stiffness 

k s =E s A s /L s (lbf/in) 

2.55xl0 4 

3.596xl0 5 

1.2xl0 6 

Critical Dynamic Pressure Ac- 

322 

140 

137 
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Now that the flutter boundaries of all five member beams have been studied, the law 
beyond them loomed. The cantilever beam with Ac- =136 is both the down and up limit 
cases of the clamped-spring stopper beam. The clamped-simply supported beam with Acr 
= 480 is the upper limit case of the connected spring stopper beam. This can be observed 
from the right-hand end of the curve in Fig. 12. Both the spring stopper beam and 
connected spring stopper beam exert constrain to the free end of the cantilever beam. So 
they will enhance the flutter boundary due to stiffening of the whole dynamic system. It's 
verified by the appearance of look-like ascending stages appeared in both Fig. 5 and Fig. 
12. But the constrain from spring stopper is unidirectional, so spring stopper beam is 
weaker than connected spring stopper beam which applies a bidirectional constrain to the 
beam. The descending of flutter boundary is then reasonable. Analysis of bar stopper 
beam proved this issue. 

In summary, the flutter boundary curve shown in Fig. 5 exposed the character of the 
original impact problem. The time domain analysis procedure was proved to be valid for 
flutter analysis of systems with impact boundary conditions. However, it seems current 
work accomplished only the prelude of a systematic work, since the practical problem has 
inevitable effects from temperature, structural damping, etc.. While the model is 
expanded to 3D plate and composite materials are considered, predictably, flow direction 
effect needs attention. 
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Fig. 2 Cantilever beam with stopper 
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II Spring Stopper ( nonlinear system) 
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IV Connected Spring Stopper 
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V Simply Supported 


Fig. 3 Various support conditions at the right end 
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Fig. 5 Flutter boundary of clamped-spring stopper beam 
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Fig. 6 Cantilever beam critical dynamic pressure 
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Time Response for k=0.4 X=160 












